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ABSTRACT 


In  the  past  several  years  there  have  been  several  large  (M^  >  7.0)  earthquakes  in  the  eastern  Mediterranean 
region  (Gulf  of  Aqaba,  Racha,  Adana,  etc.),  many  of  which  have  had  aftershock  deployments  by  local 
seismological  organizations.  In  addition  to  providing  ground  truth  data  (GT  «  5  km)  that  is  used  in 
regional  location  calibration  and  validation,  the  waveform  data  can  be  used  to  aid  in  calibrating  regional 
magnitudes,  seismic  discriminants,  and  velocity  structure.  Eor  small  regional  events  (mb  «  4.5),  a  stable, 
accurate  magnitude  is  essential  in  the  development  of  realistic  detection  threshold  curves,  proper  magnitude 
and  distance  amplitude  correction  processing,  formation  of  an  Ms:mb  discriminant,  and  accurate  yield 
determination  of  clandestine  nuclear  explosions. 

Our  approach  provides  a  stable  source  spectra  from  which  M^,  and  mb  can  be  obtained  without  regional 
magnitude  biases.  Once  calibration  corrections  are  obtained  for  earthquakes,  the  coda-derived  source 
spectra  exhibit  strong  depth-dependent  spectral  peaking  when  the  same  corrections  are  applied  to 
explosions  at  the  Nevada  Test  Site  (Mayeda  and  Walter,  1996),  chemical  explosions  in  the  recent  “Depth  of 
Burial”  experiment  in  Kazahkstan  (Myers  et  al.,  1999),  and  the  recent  nuclear  test  in  India.  Eor  events  in 
the  western  U.S.  we  found  that  total  seismic  energy,  E,  scales  as  resulting  in  more  radiated  energy 

than  would  be  expected  under  the  assumptions  of  constant  stress-drop  scaling.  Preliminary  results  for 
events  in  the  Middle  East  region  also  show  this  behavior,  which  appears  to  be  the  result  of  intermediate 
spectra  fall-off  (f  *  ^)  for  frequencies  ranging  between  ~0. 1  and  0.8  Hz  for  the  larger  events.  We  developed 
a  Seismic  Analysis  Code  (SAC)  coda  processing  command  that  reads  in  an  ASCII  flat  file  that  contains 
calibration  information  specific  for  a  station  and  surrounding  region,  then  outputs  a  coda-derived  source 
spectra,  moment  estimate,  and  energy  estimate. 

This  is  Lawrence  Livermore  National  Laboratory  document  nurrberUCRL-JC-138990. 
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OBJECTIVE 

For  sparse  local  and  regional  seismic  networks  a  stable  magnitude  is  of  utmost  importance  for  establishing 
accurate  seismicity  catalogs  and  assessing  seismic  hazard  potential.  In  the  context  of  network  capability, 
we  need  accurate  magnitudes  for  construction  of  detection  threshold  curves,  a  crucial  element  of  nuclear 
test  monitoring.  In  this  paper  we  describe  a  method  of  calibrating  a  region  for  coda  envelope  processing 
using  a  new  routine  in  SAC.  Unlike  magnitudes  such  as  Ml,  M^  and  mb  which  are  relative,  narrowband 
measures  that  often  have  regional  biases,  our  approach  can  provide  a  stable  absolute  source  spectra  that  is 
corrected  for  S-to-coda  transfer  function,  scattering,  and  anelastic  attenuation  and  site  effect.  Because  coda 
envelope  amplitude  measurements  are  not  as  sensitive  to  3-D  path  heterogeneity  and  source  radiation 
pattern,  the  resultant  spectra  can  then  be  used  for  stable  moment  estimation  (and  hence  M„),  traditional 
short -period  magnitudes  (i.e.,  nib.  Ml)  and  radiated  seismic  energy,  E.  There  have  been  numerous 
scattering  studies  over  the  past  3  decades  investigating  the  mechanisms  that  control  coda  generation  and  its 
decay  with  time  from  the  origin  and  dependence  upon  distance.  For  our  method  however,  we  only  use  an 
empirical  approach  to  fit  the  observed  envelopes  to  obtain  an  amplitude  measure  over  a  range  of  narrow 
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frequency  bands.  We  then  empirically  tie  these  amplitudes  to  independently  derived  moment  estimates 
from  calibration  events  through  a  simple  boot-strapping  procedure.  Past  studies  for  earthquakes  in  the 
western  United  States  showed  that  a  single  coda  envelope  amplitude  measurement  was  roughly  equivalent 
to  a  64  station  network  using  direct  waves  (Mayeda  &  Walter,  1996).  For  explosions  at  the  Nevada  Test 
Site  (NTS)  we  found  that  a  1  Hz  single-station  coda  measurement  was  roughly  equivalent  to  a  20  station 
network  using  direct  Lg  or  (Mayeda,  1993).  This  is  the  result  of  making  a  time  domain  measurement 
over  a  large  length  of  coda,  thereby  averaging  the  combined  effects  of  source  radiation  pattern,  path 
heterogeneity,  and  near  site  random  interference,  all  of  which  can  affect  individual  phases  such  as  Lg  and 
Pn.  The  results  of  the  frequency-dependent  calibration  have  been  used  as  input  into  a  newly  developed 
interactive  SAC  command  that  forms  narrowband  envelopes,  measures  amplitudes,  applies  region  and 
station-specific  corrections,  outputs  ASCII  source  information  (e.g.,  M„,  nib,  Aa,  E)  and  plots  the  moment- 
rate  spectra.  In  this  paper  we  describe  the  calibration  method,  then  apply  it  to  a  data  set  in  the  Middle  East 
region  to  estimate  seismic  moment  (and  M„),  total  radiated  seismic  energy,  and  dynamic  stress  drop  (e.g., 
Orowan)  (see  Eigure  1).  We  then  compare  with  previous  results  from  the  western  United  States  (Kanamori 
et  ah,  1993,  Mayeda  and  Walter,  1996). 
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Synthetic  Envelopes 

We  found  that  a  simple  analytic  form  that  resembles  the  single -scattering  model  of  Aki  (1969)  generally 
does  a  good  job  at  fitting  both  regional  and  local  coda  envelopes.  The  main  difference  is  that  we  found  that 
we  had  better  fits  by  forming  the  envelopes  relative  to  the  S  arrival,  not  the  origin  time  as  is  done  in  most 
scattering  models.  The  analytic  expression  that  we  used  to  fit  the  observed  narrowband  envelopes  at  center 
frequency /:is 


A  ,(f,  1 1  x)  =  A„  •  H(t  -  ^)  •  (t  -  •  exp(-b  •  (t  -  ^)))  (1) 

where  Ag  is  a  source  constant,  H  is  the  Heaviside  step  function,  t  is  the  time  in  seconds  from  the  origin,  x  is 
the  distance  in  kilometers,  v(v)  is  the  velocity  of  the  peak  arrival  in  km/sec  and  y  and  b  control  the  coda 
decay.  It  should  be  noted  however  that  the  more  complex  2-D  and  3-D  multiple  scattering  models  could 
also  be  used  (see  Sato  and  Eehler  review,  1999).  Eor  our  purposes  we  only  want  an  empirical  fit  to  the  data 
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over  a  range  of  distances  and  frequency  bands  and  thus  this  simple  form  is  completely  adequate  for  this 
purpose. 

Since  we  are  not  explicitly  applying  a  scattering  model  that  has  attenuation  built  into  its  formulation,  we 
need  to  find  an  empirical  distance  relation  for  the  coda  amplitudes.  This  is  preferable  since  we  do  not  want 
to  make  any  assumptions  about  homogeneous  coda  decay  in  a  region,  an  assumption  commonly  made  in 
past  coda  studies.  Instead,  we  fold  in  both  geometrical  spreading  and  attenuation  (whether  it  be  from 
scattering,  absorption,  leakage  etc.)  using  a  strictly  empirical  form.  We  reiterate  that  we  use  equation  1  as  a 
means  of  matching  the  observed  coda  envelope  shape  to  extract  a  coda  amplitude  measurement.  In 
essence,  equation  I  is  used  as  an  empirical  metric  and  thus  any  functional  relationship  could  have  been 
used  that  fits  the  data.  The  steps  described  in  this  paper  allow  one  to  measure  coda  envelopes,  which  are 
initially  in  dimensionless  units,  correct  for  distance  and  site  effects  and  tie  to  an  absolute  measure.  The 
validation  of  this  approach  is  simple.  We  verify  that  we  obtain  the  same  source  spectra  at  different  stations 
and  distances  for  the  same  event  (thus  confirming  our  empirical  distance  corrections)  and  then  verify  that 
our  inferred  moments  are  the  same  as  those  that  were  independently  determined  from  other  means  such  as 
long-period  waveform  modeling. 

Coda  Envelopes 

For  our  data  we  applied  an  8*'’  order,  zero-phase  (four  poles,  two  passes)  Butterworth  filter  for  13  narrow 
frequency  bands  ranging  between  0.03  and  8.0  Hz  for  the  two  horizontal  components.  This  range  of 
frequency  bands  was  chosen  such  that  we  could  compute  for  the  Gulf  of  Aqaba  mainshock  of 
November  22,  1995  (M„  7.1)  and  measure  the  high  frequencies  of  some  of  the  smaller  aftershocks. 

For  each  component  the  narrowband  envelope  at  center  frequency/is  of  the  form; 

e(t  lf)=V(v(t)"  +  h(t)^)  (2) 

where  v(f)  is  one  of  the  bandpass  filtered  seismograms  and  h(t)  is  its  corresponding  Hilbert  transform.  We 
then  take  the  logio  of  the  envelopes,  average  the  two  horizontal  envelopes,  then  smooth.  Fig.  (2)  shows  an 
example  of  a  broadband  trace  and  a  selection  of  envelopes  in  different  narrow  frequency  bands.  It  is 
preferable  to  use  the  horizontal  elements  since  S-waves  have  better  signal-to-noise  ratio  and  averaging  the 
two  provides  a  smoother  envelope  than  a  single  component  alone.  Of  course  this  processing  could  also  be 
done  on  a  single  component  such  as  the  vertical  but  is  less  desirable. 
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Velocity  of  Direct  Phase 

After  forming  narrowband  envelopes  for  all  the  events,  we  measure  the  (phase)  velocity  of  the  peak  S- 
arrival  for  each  frequency  band  and  plot  versus  distance.  For  the  longest  periods  the  peak  corresponds  to 
the  Rayleigh  or  Love  waves  but  because  we  empirically  correct  for  each  frequency  band,  mixing  wavetypes 
makes  no  difference  on  the  final  results  as  will  be  shown  in  the  concluding  sections.  Since  our  goal  is  to 
calibrate  both  local  and  regional  distance  earthquakes,  we  find  that  a  simple  second-order  polynomial  fits 
the  data  at  all  distances.  For  each  frequency  the  distant-dependent  velocity,  v  is  of  the  form, 

v(x  I  f)  =  velorder_0  +  velorder_l  •  x  -b  velorder  _2  •  x"  (3) 

where  velorderJO,  velorder_l  and  velorder_2  are  constants. 

Coda  Shape  Factors 

Now  that  we  know  the  velocities  for  each  frequency  band,  we  can  fit  our  simple  analytic  form  to 
find  the  coda  shape  parameter  b  which  controls  the  coda  decay  as  measured  from  the  direct  S  (or  Lg) 
arrival.  By  rearranging  equation  1  and  taking  the  logio  of  both  sides  we  get, 

logio  [Ac(f.  1 1  x)  •  (t  -  ]=  logio(A„  •  H(t  -  ))-  ^  •  (t  -  •  log(e)  (4) 


We  let  Y  vary  from  0  to  2  and  found  that  fixing  7=0.7  fit  the  data  well  for  the  entire  frequency  range. 
Coincidentally,  this  value  is  consistent  with  other  scattering  formulations  that  consider  surface  wave  and 

X 

body  wave  scattering  in  a  full-space.  By  plotting  equation  4  versus  t  ■ 


V 


V(X)^ 


,  the  slope  of  the  best 


fitting  line  is  blog(e). 

Since  we  want  to  fit  the  coda  at  all  distances  we  try  to  measure  the  coda  over  a  range  of  distances  to 
determine  if  b  is  dependent  upon  distance.  For  the  higher  frequencies  b  is  strongly  distant  dependent 
whereas  for  frequencies  below  -0.5  Flz,  the  coda  shape  factor  is  roughly  constant  and  is  smaller  (i.e., 
decays  more  slowly  than  the  higher  frequencies).  We  find  that  a  simple  linear  form  fits  the  b  values  as  a 
function  of  distance  x  spanning  both  local  and  regional  distances. 


b(x  lf)  =  b„(f)  +  r[x 


(5) 


where  bjf)  is  the  intercept  value  at  zero  distance  at  center  frequency/and  r|  is  the  slope. 

Simple  Distance  Corrections 

For  events  in  the  Dead  Sea  rift  and  Gulf  of  Aqaba  region,  the  Geophysical  Institute  of  Israel  (GII)  routinely 
performs  network  locations  and  assorted  magnitudes  such  as  mt.  Ml  and  M^.  Now  that  we  have  the 
velocities  and  coda  shape  factors  we  setA„  in  equation  (1)  to  unity,  take  the  logio  to  be  consistent  with  our 
observed  envelopes,  then  dc  shift  the  synthetic  envelopes  to  fit  the  observed  envelopes  using  an  L-1  norm. 
The  magnitude  of  the  dc  shift  is  the  non-dimensional  coda  amplitude.  This  amplitude  is  analogous  to  a 
direct  wave  measure  in  the  sense  that  a  distance  and  site  correction  still  need  to  be  made.  Measuring  the 
coda  envelope  amplitude  over  a  length  of  time  merely  provides  a  more  stable  measure  than  using  direct 
waves  which  are  of  short  duration  and  thus  considerably  more  susceptible  to  interference,  source-station 
heterogeneity,  and  source  radiation  pattern.  In  the  following  paragraph  we  empirically  find  the  distance 
dependence  for  each  frequency  band,  apply  it  to  the  data,  then  validate  by  comparing  measurements  of  the 
same  event  at  two  different  stations. 

For  each  frequency  band  we  subtract  the  GII  Ml  from  the  amplitude  and  plot  versus  distance  in  kilometers. 
We  try  to  use  events  in  a  narrow  magnitude  range  for  each  frequency  band  when  determining  the  distance 
corrections  to  avoid  potential  biases  related  to  corner  frequency  scaling.  For  example,  using  only  small 
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events  at  short  distances  and  only  large  events  at  longer  distances  would  result  in  very  biased  curve  fitting 
as  opposed  to  using  events  of  roughly  the  same  size.  Up  to  this  point  we  can  now  measure  distance- 
corrected  coda  amplitudes  from  narrowband  envelopes  but  these  are  still  in  dimensionless  units.  Figure 
(3a)  shows  a  cartoon  of  non-dimensional  coda  spectra  for  hypothetical  earthquakes  of  various  sizes  and  its 
corresponding  figure  shows  what  the  spectra  will  look  like  after  correcting  to  an  absolute  scale.  The  goal  is 
to  correct  for  S-to-coda  transfer  function  and  site  effects  to  obtain  a  moment-rate  spectra.  Again,  these 
amplitudes  are  relative  to  a  source  amplitude  of  1.0  and  are  thus  in  dimensionless  units  and  must  be 
corrected  to  an  absolute  scale.  The  following  section  describes  this  crucial  series  of  steps. 
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Figure  3a  arxl  3b  Non-dimenslonal  coda  amplitude  corrected  lor  distance  effect  only  (fig  3a)  for  a  range 
of  events  Figure  3b  shows  corrected  spectra  after  tying  to  mdependenl  moments  arw)  applying  Green's  functxx) 
corrections 


Green’s  Function  and  Seismic  Moment  Corrections 

Up  to  this  point  we  have  measured  non-dimensional  coda  envelope  amplitudes  and  corrected  for  distance. 

In  this  section  we  relate  these  values  to  an  absolute  scale,  namely  seismic  moment  to  obtain  a  moment-rate 
spectrum.  Our  only  assumption  is  that  the  S-wave  source  spectrum  is  flat  below  the  corner  frequency.  We 
purposely  underestimate  where  the  corner  frequency  lies  for  our  set  of  calibration  events  to  avoid  flattening 
our  spectra  unrealistically. 

We  used  a  1-D  reflectivity  code  to  waveform  model  a  series  of  events  in  the  region  to  estimate  seismic 
moment.  The  only  assumption  we  make  is  that  the  source  spectrum  of  these  events  is  flat  below  the  corner 
frequency.  As  shown  in  Figure  (3b),  we  add  constants  to  all  amplitudes  for  each  frequency  band  such  that 
the  seismic  moments,  in  a  least  squares  sense,  agree  with  the  waveform  modeled  results.  Since  these  are 
frequency-dependent  corrections,  independent  of  distance,  the  corrections  must  be  uniformly  applied  to  all 
distance-corrected  amplitudes.  Figure  (3b)  shows  that  with  conservative  estimates  of  the  corner  frequency, 
we  can  flatten  the  spectra  for  a  range  of  event  sizes.  For  events  that  are  too  small  to  be  waveform  modeled 
(<  ~Mw  3.5),  we  can  stiU  estimate  moment  using  periods  less  than  a  few  Hertz.  For  especially  small  events 
(~Mw  2.5)  we  can  assume  they  have  flat  source  spectrum  at  least  to  10  Hz.  This  use  of  the  very  small 
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events  is  essentially  a  Green’s  function  correction  that  allows  us  to  extend  the  stable  coda-derived  source 
spectrum  to  much  smaller  events.  Therefore,  this  allows  us  to  find  the  amplitude  corrections  for  the 
frequencies  between  ~2.0  and  8  Hz.  Though  in  figures  (3a  &  3b)  we  use  just  a  few  events  for  illustrative 
purposes,  for  the  actual  data  we  used  a  number  of  events  of  the  same  relative  size  to  come  up  with  the  best 
set  of  moment  corrections.  In  our  case  we  instrument  deconvolved  the  traces  beforehand  and  converted  to 
velocity.  However,  if  the  instrument  did  not  change  over  the  range  of  time  of  study,  the  moment 
corrections  (which  include  site  and  S-to-coda  transfer  function)  could  have  also  included  the  instrument 
response. 


Aqaba  Spectra 


Frequency  (Hz) 


We  have  good  confidence  in  our  magnitudes  for  events  that  are  too 
smaii  to  be  waveform  modeied. 

Figure  4  Selected  moment-fate  spectra  lor  the  Quit  of  Aoaba  sequence 


TESTING  &  VALIDATION 

Figure  (4)  shows  selected  source  spectra  for  events  in  the  Gulf  of  Aqaba  region  along  with  some  events 
from  the  Dead  Sea  rift  at  station  BGIO.  We  test  for  consistency  by  comparing  spectra  of  the  same  event  at 
two  different  stations,  in  our  case,  stations  BGIO  in  Israel  and  KEG  in  Egypt.  Eigure  (5)  shows  spectra  for 
two  events,  one  is  nodal  at  BGIO  and  the  other  is  not.  In  both  cases,  the  spectra  at  KEG  and  BGIO  are 
virtually  identical,  despite  a  significant  source  radiation  pattern  difference.  If  we  compare  individual 
amplitudes  for  a  large  number  of  events,  we  also  see  that  the  inter-station  variation  is  also  very  small, 
comparable  to  what  was  observed  at  NTS. 
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•  Coda  ts  not  Mnattive  to  aourco 
radiation  pattern. 

•  Coda-darived  magnitudes  are  possibte 
for  significantly  smaller  events. 

•  Results  are  comparable  to  the 
Western  United  States. 

•  ¥w  Is  tied  to  a  physical  measure 
of  earthquake  size  and  is  unbiased. 


Figure  5  Spectm  recorded  at  KEG  and  BGIO  kx  the  same  two  events  show  tMe  vanstion 
Note  The  November  23. 1905  attefshockw  nodal  at  BGlO  but  the  coda-derrved 
spectra  Is  the  same  as  the  spectra  at  KEG 


Finally,  we  verify  that  our  absolute  amplitudes  are  correct  by  computing  M*  from  the  lowest  frequencies  of 
our  measured  spectra  and  comparing  against  independent  M^’s.  Figure  (6)  shows  that  there  is  good 
agreement  between  the  two  approaches. 
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agreed  with  wavetomi  modMmg  results 
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CONCLUSIONS  &  RECOMMENDATIONS 


Now  that  we  have  verified  that  our  spectra  are  completely  corrected  back  to  the  S-wave  source,  we  can 
begin  to  address  a  long-standing  problem  in  observational  seismology,  namely  accurate  estimation  of  the 
radiated  seismic  energy.  For  events  for  which  we  observe  a  flat,  long-period  level  (proportional  to 
moment)  and  decaying  high  frequencies,  we  need  to  extrapolate  our  observed  spectra  to  f=0  and  f=  oo  . 
Though  a  majority  of  the  energy  is  near  the  corner  frequency,  for  completeness  we  must  do  this 
extrapolation.  For  the  long  periods  we  extend  our  spectra  to  f=0  Hz  and  for  the  high  frequencies,  we 
assume  an  omega  square  decay  to  f=  °o  .  Next,  we  convert  our  spectra  to  velocity  by  dividing  by  omega, 
then  square  and  integrate.  Only  those  events  for  which  we  have  measured  at  least  70%  of  the  total  S-wave 
energy  are  used.  Finally,  we  assume  a  P-wave  contribution  of  energy  by  multiplying  by  1.07.  From  Figure 
7  we  see  that  the  Gulf  of  Aqaba  events  follow  the  same  trend  as  those  from  the  western  United  States. 

We  have  successfully  transported  a  methodology  for  calibrating  stations  for  regional  magnitude  to 
International  Monitoring  System  (IMS)  stations  in  the  Middle  East  region.  Resultant  source  spectra  show 
very  small  interstation  variation  as  was  found  in  the  western  U.S.  and  will  allow  accurate  estimates  of  Mw 
and  mb  for  small  regional  events  that  cannot  be  measured  teleseismically.  We  are  systematically 
calibrating  all  the  IMS  stations  in  the  region  and  are  tying  regional  mb’s  to  the  teleseismic  mb.  The  coda- 
derived  Mw  estimates  are  being  used  as  input  to  the  joint  LANL/LLNL  MDAC  code  for  discrimination 
calibration.  The  coda  methodology  has  been  incorporated  into  the  newest  version  of  SAC  that  will  allow 
amplitude  measurements  on  narrowband  coda  envelopes  given  a  calibration  flatfile. 

Radiated  Seismic  Energy  and  Orowan  Stress  Drop 


Figure  7.  Energy-moment  relationship  tor  Aqaba  events  is  in  good  agreement  with 
results  from  the  western  United  States. 
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